package main
import (
"os"
"fmt"
"strconv"
"math"
"rand"
)

func Rand() float64{
//    return 0.4
    return rand.Float64()
}

type Vec struct {
    x,y,z float64
}

func (a *Vec) String() string{
    return fmt.Sprintf("[%f,%f,%f]",a.x,a.y,a.z)
}

func (a Vec) add(b Vec) Vec{
    return Vec{a.x+b.x,a.y+b.y,a.z+b.z}
}

func (a Vec) dec(b Vec) Vec{
    return Vec{a.x-b.x,a.y-b.y,a.z-b.z}
}

func (a Vec) mul (b float64) Vec{
    return Vec{a.x*b,a.y*b,a.z*b}
}

func (a Vec) mulv(b Vec) Vec{
    return Vec{a.x*b.x,a.y*b.y,a.z*b.z}
}

func (a *Vec)norm() Vec{
    *a =  (*a).mul(1/math.Sqrt(a.x*a.x + a.y*a.y + a.z*a.z))
    return *a
}

func (a Vec)dot(b Vec) float64{
    return a.x*b.x + a.y*b.y + a.z*b.z
}

func (a Vec)mod(b Vec) Vec{
    return Vec{a.y*b.z - a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x}
}

type Ray struct{
    o,d Vec
}

func (r *Ray) String() string{
     return fmt.Sprintf("{%s,%s}",r.o,r.d)
}

type Refl_t int
const (
    DIFF Refl_t = iota
    SPEC
    REFR
)

type Sphere struct {
    rad float64
    p,e,c Vec
    refl Refl_t
}

func (s *Sphere) intersect (r *Ray) float64 {
    op := s.p.dec(r.o)
    var eps float64
    eps = 1e-4
    b := op.dot(r.d)
    det := b*b - op.dot(op) + s.rad *s.rad
    if det < 0{
        return 0
    }
    det = math.Sqrt(det)
    t := b-det
    if t > eps {
        return t
    }else  {
        t = b + det
        if t > eps {
            return t
        }
    }
    return 0
}

func clamp (x float64) float64{
    if x < 0 { return 0}
    if x > 1 { return 1}
    return x
}

func toInt (x float64) int {
    return int(math.Pow(clamp(x),1/2.2)*255 + .5)
}

type  SphereArrary []Sphere
func (spheres *SphereArrary) intersect (r *Ray,t *float64,id *int) bool{
    var d,inf float64
    *t = 1e20
    inf  = 1e20
    for i:=len(*spheres)-1;i>=0;i--{
        d = (*spheres)[i].intersect(r)
        if d > 0 && d < *t{
            *t = d
            *id = i
        }
    }
    return *t < inf
}

func (spheres *SphereArrary)radiance (r *Ray,depth int,Xi []uint16) Vec{
    var t float64
    id := 0
    if(!spheres.intersect(r,&t,&id)) {
        return Vec{}
    }
    
    obj := &(*spheres)[id]
    var x,nl,f Vec
    var n Vec
    x = r.o.add(r.d.mul(t))
    n = x.dec(obj.p)
    n = (&n).norm()
    if n.dot(r.d) < 0{
        nl = n
    }else{
        nl = n.mul(-1)
    }
    f = obj.c
    //max refl
    var p float64
    if f.x>f.y && f.x > f.z {
        p = f.x
    }else if f.y > f.z {
        p = f.y
    }else {
        p = f.z
    }

    depth++
    if depth > 5{
        if Rand() < p {
            f = f.mul(1/p)
        }else{
            return obj.e
        }
    }
    if obj.refl == DIFF {
        r1 := 2*math.Pi*Rand()
        r2 := Rand()
        r2s := math.Sqrt(r2)
        var w,u,v Vec
        w = nl
        if  math.Fabs(w.x) > .1 {
            u = Vec{0,1,0}
        }else{
            u = (Vec{1,0,0})
        }
	u = u.mod(w)
	u = (&u).norm()
        v = w.mod(u)
        d := (u.mul(math.Cos(r1)).mul(r2s).add(v.mul(math.Sin(r1)).mul(r2s)).add(w.mul(math.Sqrt(1-r2))))
	d = (&d).norm()

        return obj.e.add(f.mulv(spheres.radiance(&Ray{x,d},depth,Xi)))
    }else if (obj.refl == SPEC){
        return obj.e.add(f.mulv(spheres.radiance(&Ray{x,r.d.dec(n.mul(2).mul(n.dot(r.d)))},depth,Xi)))
    }
    reflRay := &Ray{x,r.d.dec( n.mul(2).mul(n.dot(r.d)))};
    into := n.dot(nl) > 0;
    var nc,nt,nnt,ddn,cos2t float64
    nc = 1
    nt = 1.5
    if into {
        nnt = nc/nt
    }else{
        nnt = nt/nc
    }
    ddn =  r.d.dot(nl)
    cos2t = 1-nnt*nnt*(1-ddn*ddn)
    if cos2t < 0{
        return obj.e.add(f.mulv(spheres.radiance(reflRay,depth,Xi)))
    }

    var times float64
    times = 1
    if !into {
        times = -1
    }

    tdir := (r.d.mul(nnt).dec(n.mul((times * (ddn*nnt + math.Sqrt(cos2t))))))
    tdir = (&tdir).norm()
    a := nt - nc
    b := nt + nc
    R0 := a*a/(b*b)
    var c float64
    if into {
        c = -ddn
    }else{
        c= tdir.dot(n)
    }
    c = 1 - c
    Re := R0 + (1-R0)*c*c*c*c*c
    Tr := 1 - Re
    P :=.25 + .5*Re
    RP := Re/P
    TP := Tr/(1-P)
    var tmp_vec Vec
    if depth > 2 {
        if(Rand() < P ){
            tmp_vec = spheres.radiance(reflRay,depth,Xi).mul(RP)
        }else{
            tmp_vec = spheres.radiance(&Ray{x,tdir},depth,Xi).mul(TP)
        }
    }else{
        tmp_vec = spheres.radiance(reflRay,depth,Xi).mul(Re).add(spheres.radiance(&Ray{x,tdir},depth,Xi).mul(Tr))
    }
    return obj.e.add( f.mulv(tmp_vec))
}



func main(){
    spheres  := SphereArrary{
        //Scene: radius, position, emission, color, material
        Sphere{1e5, Vec{ 1e5+1,40.8,81.6}, Vec{},Vec{.75,.25,.25},DIFF},//Left
        Sphere{1e5, Vec{-1e5+99,40.8,81.6},Vec{},Vec{.25,.25,.75},DIFF},//Rght
        Sphere{1e5, Vec{50,40.8, 1e5},     Vec{},Vec{.75,.75,.75},DIFF},//Back
        Sphere{1e5, Vec{50,40.8,-1e5+170}, Vec{},Vec{},           DIFF},//Frnt
        Sphere{1e5, Vec{50, 1e5, 81.6},    Vec{},Vec{.75,.75,.75},DIFF},//Botm
        Sphere{1e5, Vec{50,-1e5+81.6,81.6},Vec{},Vec{.75,.75,.75},DIFF},//Top
        Sphere{16.5,Vec{27,16.5,47},       Vec{},(Vec{1,1,1}).mul(.999), SPEC},//Mirr
        Sphere{16.5,Vec{73,16.5,78},       Vec{},(Vec{1,1,1}).mul(.999), REFR},//Glas
        Sphere{600, Vec{50,681.6-.27,81.6},Vec{12,12,12},  Vec{}, DIFF}, //Lite
    }


    w := 1024
    h := 768
    samps := 1
    if len(os.Args) == 2{
        samps,_ = strconv.Atoi(os.Args[1])
        samps = samps/4
    }
    cam := &Ray{Vec{50,52,295.6},(&Vec{0,-0.042612,-1}).norm()}
    cx := Vec{float64(w)*.5135/float64(h),0,0}
    cy := cx.mod(cam.d) 
    cy = (&cy).norm().mul(.5135)
    var r Vec
    c := make([]Vec,w*h)
    for y:=0;y<h;y++ {
        fmt.Printf("\rRendering (%d spp) %5.2f%%",samps*4,float64(100)*float64(y)/float64(h-1));
        for x:=0;x < w;x++ {
            sy:= 0
            for i := (h-y-1)*w +x ; sy < 2;sy++{
                var Xi [3]uint16
                for sx:=0;sx<2;sx++{
                    for s:=0;s<samps;s++{
                        r=Vec{}
                        r1 := 2*Rand()
                        var dx,dy float64
                        if r1 < 1{
                            dx = math.Sqrt(r1)-1
                        }else{
                            dx = 1- math.Sqrt(2-r1)
                        }
                        r2 := 2*Rand()
                        if r2 < 1{
                            dy = math.Sqrt(r2)-1
                        }else{
                            dy = 1-math.Sqrt(2-r2)
                        }

                        d := cx.mul((((float64(sx)+.5+dx)/2+float64(x))/float64(w) -.5)).add(
                            cy.mul(((float64(sy)+.5+dy)/2+float64(y))/float64(h) -.5) ).add(cam.d)
                        d.norm()
                        r = r.add(spheres.radiance(&Ray{cam.o.add(d.mul(140)),d},0,&Xi).mul(1./float64(samps)))
                    }
                    c[i] = c[i].add((Vec{clamp(r.x),clamp(r.y),clamp(r.z)}).mul(.25))
                }
            }
        }
    }
    f,_ := os.Open("image.ppm",os.O_WRONLY | os.O_CREATE,0600)
    defer f.Close()
    f.WriteString(fmt.Sprintf("P3\n%d %d\n%d\n",w,h,255))
    for i:=0;i<w*h;i++{
        f.WriteString(fmt.Sprintf("%d %d %d ",toInt(c[i].x),toInt(c[i].y),toInt(c[i].z)))
    }
    return
}
